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ABSTRACT 



Techniques for the statistical analysis of the Lyman a forest in high redshift 
quasars are developed, and applied to the low resolution (25 A) spectra of 29 of 
the 33 quasars in the Schneider-Schmidt-Gunn (SSG) sample. We extrapolate 
Q \ each quasar's continuum shortward of Lyman a emission, then consider each 

spectral bin of each quasar to be an (approximately) independent measurement 
fT) \ of the absorption due to the Lyman a clouds. With several thousand such 



measurements thus available, we can obtain good determinations of some 
interesting properties of clouds in the redshift range 2.5 < z < 4.3 without 
actually resolving any single cloud. We find that the mean absorption increases 
with z approximately as a power law (1 + z) 1+1 with 7 = 2.46 ± 0.37. The mean 



ratio of Lyman a to Lyman j3 absorption in the clouds is 0.476 ± 0.054. We 
also detect, and obtain ratios, for Lyman (3, 7, and possibly e. We are also able 
to quantify the fluctuations of the absorption around its mean, and find that 
these are comparable to, or perhaps slightly larger than, that expected from an 
uncorrelated distribution of clouds. The techniques in this paper, which include 
the use of bootstrap resampling of the quasar sample to obtain estimated errors 
and error covariances, and a mathematical treatment of absorption from a 
(possibly non-uniform) stochastic distribution of lines, should be applicable to 
future, more extensive, data sets. 



Subject headings: cosmology: observations - quasars - intergalactic medium 
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1. Introduction 

Excepting only the quasars themselves, no observable objects have more potential 
for revealing the quantitative nature of the Universe at early times than do the Lyman a 
clouds, seen in absorption against the UV continuum of background quasars. Indeed, one 
might argue that the importance of the Lyman a clouds, at redshifts z > 2.5, in some 
respects exceeds that of the quasars, since the high-redshift quasars are "unusual" objects, 
presumably associated with the extreme statistical tail of structure formation, while the 
Lyman a clouds are probably more "typical" representatives of the state of baryonic matter, 
at least in the redshift range 2.5 < z < 4.3. 

From the time of their first discovery (Lynds 1971) and first detailed analysis (Sargent 
et al. 1980), the Lyman a clouds have principally been studied by the spectroscopic 
identification, and counting, of individual clouds along a quasar line of sight. (For a review, 
see, e.g., Sargent 1988a.) Since the cloud absorption lines are narrow, with equivalent 
widths on the order of 1 A or less, observations with moderate-to-high spectral resolving 
power are required. Indeed, some exceptionally high quality QSO spectra (resolution 
^ 0.1 A) have been obtained for use in direct studies of the cloud line profiles and velocity 
dispersion parameters (Pettini et al. 1990, Carswell et al. 1991). 

At the opposite extreme lies the set of techniques that study the gross relative 
depression of the quasar continuum spectrum shortward of Lyman a (that is, emitted 
wavelength X em < 1216 A) by the aggregate effect of many clouds along the line of sight at 
redshifts smaller than the emission redshift z em (Oke and Korycansky 1982; Bechtold et 
al. 1984; Schneider et al. 1989a,b, 1991; Giallongo and Cristiani 1990; Jenkins and Ostriker, 
1991). Typical of these techniques is the association of a single number Da with each 
quasar line of sight, the fractional absorption averaged over a broad band safely between 
the Lyman a and Lyman (3 emission wavelengths, say 1050A < A em < 1170A. 

Despite the obvious fundamental importance of understanding the nature of the 
intergalactic medium at high redshifts (z > 3.5, say), there have been very few high-redshift 
studies (Schneider et al. 1989b, 1991; Jenkins and Ostriker, 1991), and these have all used 
the Da approach - because the quasars available for study are so faint. 

In this paper we develop and apply a new technique, intermediate between the above 
extremes (though in most ways closer to the latter). Our data set consists of the carefully 
calibrated, low resolution spectra of 33 high-redshift quasars obtained by Schneider, 
Schmidt, and Gunn (1991; hereafter "SSG"). This sample provides the best available 
data set for z > 3.5 because (i) it contains all of the quasars with published redshifts 
z > 3.85 except PC 1247+3406 (z = 4.897, too large for this study), (ii) the spectra were 
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all obtained with the same instrumentation and thus have similar resolution and noise 
properties, and (hi) the sample contains some redshifts below z = 3.5, so there is some 
overlap with "detailed" studies. Our new approach (but cf. Webb et al. 1992) is to consider 
each separately resolved wavelength of each individual quasar to be an independent, albeit 
noisy, measurement of the cloud absorption in a particular three-dimensional volume of the 
Universe. (For A em < 1025A, the emission wavelength of Lyman (3, each volume consists of 
more than one disjoint piece.) 

We will find that, for the 25A spectral resolution of the SSG data, there are typically 
several or dozens, but not hundreds, of significantly absorbing clouds in each resolution 
element. Since, however, we have not 29 measured numbers, but rather several thousand, it 
is possible to untangle the statistics of the overlapping clouds to quite a remarkable degree. 
We will see, for example, that it is possible to obtain meaningful measurements of not only 
Lyman a absorption, but also Ly-/3 (1025A), Ly-7 (972 A), Ly-5 (949 A), and possibly Ly-e 
(937A). 

Jenkins and Ostriker (1991) showed that useful information could potentially be 
obtained from the detailed distribution of absorption values seen (at low or moderate 
resolution) in the Lyman a forest of a single quasar's line of sight. We extend that important 
idea in this paper, and show that, from statistics of the forest's fluctuating absorption, one 
can derive fairly precise measurements of several statistical quantities associated with the 
Lyman a cloud distribution. 

The output of this paper is a set of techniques for measuring several quantitative 
statistical properties of the Lyman a cloud distribution, along with detailed uncertainty 
estimates (including cross-correlations), and also the results obtained by applying these 
techniques to the SSG sample quasars at redshifts 2.6 < z a } }S < 4.2. Paper II of this 
series will show that these measured values are already precise enough to impose strong 
constraints on the physical nature of individual clouds and on the distribution function of 
the cloud population. Using these derived observational constraints, Paper III will be able 
to confront the grander cosmological questions associated with the clouds: their origin, 
confinement mechanism, and implications for the formation of galaxies and large scale 
structure. 

The focus of this paper, and later papers in this series, is on those Lyman a clouds 
at sufficiently high redshift to be considered (at least potentially) primordial cosmological 
objects, not associated with galaxies and minimally polluted by stellar element production. 
It is already clear (see, e.g, Sargent 1988b, Bahcall et al. 1992) that the detailed properties 
of quasar absorption clouds at low redshifts may differ significantly from those at high 
redshifts. We do not expect the results of this paper to apply at low redshifts; the question 
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"how low is low?" , obviously an interesting one, will be considered in a later paper. 

In this paper, §2 discusses how we fit the individual quasar spectra to obtain the 
underlying continua. In §3 we estimate the mean optical depth in Lyman a absorption 
as a function of redshift z. §4 extends this estimate to the other accessible Lyman lines. 
§5 analyzes the fluctuation statistics of the absorption measurements. §6 summarizes our 
conclusions. 

In Appendix A we analyze in detail a class of statistical models for the distribution of 
spectral lines in the Lyman forest. These are similar to the usual Poisson processes, but 
allow for a nonuniform underlying density of points. We obtain results for the mean of the 
transmission coefficient and also for its correlation function. In Appendix B the results of 
Appendix A are specialized for the particular applications of this paper. 



2. Estimating the Continuum 

The input data to this study are the calibrated SSG spectra of 33 quasars in the 
redshift range of 3.1 < z ern < 4.8. The data are in the form of corrected (for atmospheric 
absorption, reddening, and instrumental response) fluxes, in 10 A bins, from 4310 A to 9500 
A observed wavelength. The resolution of the measurements is about 25 A, so neighboring 
bins are not independent. (We make use of this oversampling below.) Details of the 
observations and data reduction procedures are described in SSG. 

Three of the SSG quasars are broad absorption line (BAL) quasars. We eliminate 
these from the sample. At a later stage of our processing, we also eliminate one additional 
quasar, the highest redshift member of the sample, because its available data do not 
adequately determine certain fitted parameters (see below). The 29 remaining quasars, 
along with their redshifts and monochromatic AB magnitudes at emitted wavelength 1450 
A as determined by SSG, are listed in Table 1. (Emitted wavelength 1450 A, in the gap 
between the SilV/OIV blend at 1400 A and CIV at 1549 A, is chosen as a relatively clean 
measurement of the underlying quasar continuum.) 

While quasar studies frequently focus on differences among quasar spectra, for our 
purposes it is equally important to take note of the similarity of all the spectra in this 
sample. Figure 1 plots all of the SSG measurements of all 29 quasars, as a function of 
emitted wavelength in the range 930 A to 2200 A, normalized to an (arbitrary) common 
AB U50 magnitude of 18. One sees clearly all of the spectral features present in low redshift 
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Quasar 


z 


^4-Bi45o 


Quasar 


z 




PC2344+0124 


3.143 


19.10 


PC2047+0123 


3.799 


19.23 


PC0056+0125 


3.149 


18.41 


PC1643+4631B 


3.831 


20.35 


PC1601+3754 


3.188 


19.68 


PC1301+4747 


4.004 


21.32 


PC2132+0126 


3.194 


19.69 


Q0046-293 


4.014 


19.26 


PC1605+4631 


3.203 


19.79 


PC0910+5625 


4.035 


20.86 


PC0118+0119 


3.241 


19.19 


Q0101-304 


4.072 


19.98 


PC2226+0216 


3.273 


18.80 


PC2331+0216 


4.093 


19.84 


PC0234+0120 


3.300 


19.96 


Q0000-26 


4.098 


17.46 


PC0344+0222 


3.377 


20.09 


PC0104+0215 


4.171 


19.67 


PC1619+4631 


3.471 


20.54 


PC0751+5623 


4.281 


19.65 


PC1548+4637 


3.544 


19.07 


PC0307+0222 


4.379 


19.92 


PC0345+0130 


3.638 


19.49 


Q2203+29 


4.399 


20.41 


PC1640+4628 


3.700 


19.29 


PC1233+4752 


4.447 


20.11 


PC1643+4631A 


3.790 


20.05 


PC0953+4749 


4.457 


19.09 


PC0131+0120 


3.792 


19.08 









Table 1: The 29 quasars studied in this paper. Redshifts and corrected continuum AB 
magnitudes at emitted wavelength 1450 A are from Schneider, Schmidt, and Gunn (1991). 
Four SSG quasars, the three broad absorption line quasars, and quasar PCI 158+4635 at 
redshift z = 4.733, are omitted from this study (see text). 
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quasars (see, e.g., Francis, et al. 1991). The Lyman a depression shortward of 1216 A is 
seen with equal clarity. 

While the logarithmic slope of the underlying continuum indisputably varies from 
quasar to quasar, it is striking, in Figure 1, that this is a relatively small effect over 
the wavelength range shown. The reliability of our results will depend somewhat on the 
accuracy with which we are able to extrapolate each quasar's underlying continuum down 
to ~ 930A. Fi gure 1 shows that this extrapolation, amounting to less than 25% of the range 
of abscissa shown, is not so daunting as one might think. 

To fit for the continua, we use all the available measurements between A em = 1250A 
and A em = 2200A. For each quasar in the sample, we proceed as follows. First, we estimate 
a relative statistical error associated with each measured data point. This is done by (a) 
calculating the variance of each point with its immediately adjacent neighbors, and (b) 
convolving this series of raw variances with a Gaussian profile of width about 12 bins 
(FWHM). (Here is where we exploit the fact that the data is oversampled.) It is not 
important that this estimate of statistical error be correctly normalized, but only that it 
reflect, in a general way, the relative weight to be assigned to individual measurements in 
the next stage of fitting. 

Now, we fit the data to the following 22-parameter linear model in the emitted rest 
frame, 

10 

F(A) = C 1/2 A 1/2 + C 1 A + £ 

i=l 

for the parameters C1/2, C±, A\ . . . A\§ and B\. . . Bio. Parameters C1/2 and C\ characterize 
the underlying continuum with two degrees of freedom; over our limited wavelength range 
they contain information equivalent to a magnitude and a spectral index. The constant 
values Aj and Wi, i = 1, . . . , 10, are the wavelengths and nominal widths of the 10 fitted 
lines, and are given in Table 2. Parameters A\ . . . A\$ are the fitted strengths of the emission 
lines, while B\. . . B w parametrize the ratio of the fitted width to the nominal width. The 
nominal widths in Table 2 were obtained by fitting a nonlinear model to the composite 
spectrum of Figure 1. 

The functional form of equation ([J) is perhaps slightly unconventional, and is 
motivated by the desire that it be linear in all parameters. More conventionally, one might 
fit nonlinearly for a continuum magnitude and spectral index and for the intensities and 
widths of Gaussian-profile lines. However, with the present noisy data, and the virtual 
certainty that there are other features in the data besides those modeled (especially 
longward of 1700 A), we found such nonlinear fits to be quite fussy and to require 
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considerable user intervention. Since we are not here interested in details of the shapes of 
lines in their wings, but we are interested in treating the entire data set in as statistically 
homogeneous manner as possible, the relative robustness of a purely linear fit is desirable. 

It is at this stage that we eliminate quasar PC1158+4635 at redshift z = 4.733 from 
the sample: Because of the large redshift, its longer wavelength emission features are not 
well measured, and the above fit is ill determined. While we could easily fit its continuum 
with a smaller number of parameters, we instead choose to preserve the homogeneity of 
processing procedure by eliminating the quasar entirely. 

Figure 2 shows all the data in the rest frame range 1250 A to 1800 A, and the fitted 
model spectra, individually for the 29 quasars. Also plotted for each quasar in the figure is 
the adopted underlying continuum. In many cases this is simply the function 

F(X) = C 1/2 A 1/2 + dA (2) 

where the C's are fitted parameters. In other cases, however, while the full model of 
equation ([TJ) is quite well determined, its dissection into equation (0) plus a remainder is 
quite degenerate numerically, as evidenced by a nearly degenerate family of fits that trade 
off variations in C\ii and C\ against unphysical (e.g., negative) values for the line strengths. 
In these cases we refit (essentially by eye, although the process could be automated) a 
continuum of the form of equation (Q) to the output fitted curve of ([]]), requiring the 
continuum to lie below the model at selected wavelengths between emission lines. In 
virtually all cases this refitting is quite well determined, suggesting that one could readily 
replace our procedure by a single linear fit with positivity constraints; we have not, however, 
implemented this. 

The relation between equation (Q) and the index a of the more familiar power-law 
continuum model, 

fu oc u a (3) 

at some fiducial wavelength Ao is 

0.5C 1/2 + CX /2 m 
Cl/2 + CiA 

We find a good correlation (e.g., at Ao = 1450 A) between the spectral indices thus obtained 
and those reported in SSG; however, our a values are systematically larger (smaller in 
magnitude) by a few tenths, with a mean in the range —0.5 to —0.6 reported by Richstone 
and Schmidt (1980), Steidel and Sargent (1987), Warren et al. (1991), and others. 

Adopting these continuum models, and extrapolating them from 1250 A down to 930 
A, we obtain the results shown, individually by quasar, in Figure 3. That figure shows all 
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SSG measurements shortward of emitted wavelength 1200 A in relation to their respective 
extrapolated continua. All of the rest of this paper is derived from these data. No further 
use is made of the data longward of 1200 A. 



3. Mean Optical Depth as a Function of Redshift 

The emitted wavelength range 1050 A to 1170 A is substantially uncontaminated 
by either Lyman a or Lyman (3 emission, and also (since it is longward of Lyman (3) 
uncontaminated by Lyman (3 absorption by intervening clouds. While this range could in 
principle be contaminated by C IV and Mg II absorption (see, e.g., Meyer and York 1987), 
or by quasar emission lines in this wavelength region (but see below for evidence against 
this), we will adopt the conventional assumption that these effects are likely negligible 
at high redshifts; probably one is seeing a virtually pure sample of Lyman a absorption 
by clouds. The data in this range from the 29 SSG quasars are plotted in Figure 4. 
The ordinate is the ratio of observed flux to extrapolated continuum, i.e., the fractional 
transmitted flux, or transmission. The abscissa of the left panel is emitted wavelength. One 
sees a widely scattered distribution of points. The abscissa of the right panel is observed 
wavelength or (equivalently, see top scale) absorption redshift. Here one sees a strong trend 
with redshift. Note that the right panel plots only data points which satisfy the wavelength 
cuts of the left panel, not the vastly larger number of individual SSG measurements with 
observed wavelength in the range shown (4300 A to 6500 A). 

There are in fact a few data points (not shown) with transmission greater than 1, i.e., 
with observed fluxes above the extrapolated continuum, and a few points less than zero, i.e., 
where SSG's background subtraction gives negative results. These are obviously defective 
values; to eliminate some less obvious, but likely defective, values, we adopt slightly tighter 
cuts on the data, and accept points with transmission between 0.1 and 0.9. (We have 
verified that our results are not sensitive to the values of these cutoffs.) There are 1596 
surviving points in the sample. 

Converting transmission to an equivalent optical depth (by taking the negative 
logarithm), we fit the data shown in Figure 4 a model of the conventional form 

T a (z) = A(l + zy +1 (5) 

and obtain the best-fit values A = 0.0037 and 7 = 2.46. The reason for defining the 
exponent to be 7 + 1 is so that our 7 is directly comparable to its conventional usage as 
the exponent in the number distribution of clouds (see equations 2.1-2.2 of Jenkins and 
Ostriker 1991). 
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3.1. The Bootstrap Resampling Method 



As a digression, we here need to discuss in some detail how we obtain error estimates 
(or variances) on the quantities A and 7, since we will follow a similar paradigm in obtaining 
variances and covariances for various further quantities in later sections. There are various 
reasons why the formal errors that come out of the fitting procedure are meaningless: We 
don't have good error estimates on the individual measurements. (The error estimates 
used in fitting the continuum were relative, not absolute.) The individual measurements 
shown in Figure 4 are not statistically independent, both because there is more than 
one measurement per resolution bin (adjacent measurements sample overlapping cloud 
populations), and because the points associated with a single quasar share a common source 
of error in the determination of that quasar's extrapolated continuum. 

A seemingly simple, yet very powerful, method for estimating the errors is by the 
statistical bootstrap (or resampling) method (Efron 1982, Efron and Tibshirani 1986, Press 
et al. 1992). From the sample of 29 quasars, we generated repeated resampled sets of 29 
quasars by drawing randomly with replacement. The typical set will thus have on the order 
of 12 duplicated quasars. No matter; we repeat exactly the data reduction which produced 
the original determinations of A and 7. The population of resulting values for A and 7 can 
be shown (with certain technical assumptions which need not concern us) to be distributed 
around the original determinations in the same way, both variance and covariance, that the 
original determinations ought to be distributed with respect to the true value. 

More generally, and to be applied below, suppose that we have some statistical 
procedure that determines the M quantities Ri, % — 1, . . . , M, and suppose that we make 
N resamplings. We denote the kth determination of the ith quantity by Rf , k — 0, . . . , N, 
where k = uses the full sample and k = 1, . . . , N are the resamplings. Then the best 
estimates for the R^s are given by R®, i = 1, . . . , M. The standard errors <r(Ri) are 
estimated by 




k=i 



(6) 



where 



Ri 




(7) 



The correlation matrix among the R® 



s is estimated by 




k=i 



(8) 
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An important quantity is the inverse correlation matrix , since a set of model predictions 
R* for the quantities can be compared to the measured values R® by calculating 

M 

X 2 = E^-^^-i?) (9) 

which should be chi-square distributed with M degrees of freedom (see, e.g., Rybicki and 
Press 1992). Another sometimes useful set of quantities are the coefficients of correlation 
among the individual measurements, given by = Cij/faaj). 

We can now return to our particular example: To obtain error bars on our measurement 
of 7 and A, we performed 100 resamplings, obtaining the results 

7 = 2.46 ±0.37 A = 0.0037 ± 0.0024 (10) 

where the errors are 1-er. The parameters A and 7 turn out to be very highly correlated. If 
we consider 7 as the more fundamental quantity, then almost all of the error in A can be 
moved to the determination of 7, as 

A = 0.0175 - 0.00567 ± 0.0002 (11) 

The solid line in Figure 4 shows the mean absorption of equations (pD-flTO"), transformed 



back to an ordinate of transmission. The shaded band shows the result of changing 7 by 
±lcr, with A simultaneously changed by equation (|i~ip . 



The values of A and 7 that we obtain are consistent with the results of Jenkins and 
Ostriker (1991), who measure an aggregate continuum depression for each quasar 
sampled, and in good agreement with previous determinations which relied on the counting 
of individual lines. For example, Murdoch et al. (1986) obtained 7 = 2.31 ± 0.40; Bajtlik et 
al. (1988) obtained 2.36 ±0.40; Lu, Wolfe, and Turnshek (1991) obtained 2.37 ±0.26. What 
is new in this investigation is the extension to higher redshifts (Table 1), and the use of a 
large quantity of low resolution data without line counting and without aggregation into 
broad spectral bands. Note that the large values of 7 at high redshift are quite different 
from those observed at low redshift (e.g., Bahcall et al. 1993). In §5, below, and in Paper 
II, we will see the additional information that comes from an absolute determination of the 
constant A, and from our ability (not present in aggregate determinations of Da) to look 
at the distribution of individual points around the fitted mean values. 



4. Relative Absorption Strengths of Lyman /3, 7, 5, and e 
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In Figure 5 we replot all the data shortward of X em = 1200A, but now removing 
(individually for each datum at its own appropriate absorption redshift) the effect of the 
mean Lyman a absorption, equations (W) and flTU|). The ordinate is the same as Figure 



4, namely ratio of observed flux to the extrapolated continuum. The shaded band in the 
figure is simply a moving-window average; its fluctuation scale is simply an artifact of the 
chosen window size. 

One sees that, shortward of the remaining tail of Lyman a emission and down to 1025 
A, the corrected fluxes indeed scatter around flux ratio unity. There is clear O VI emission 
around 1032A, with perhaps some contribution from Lyman (3, then a visible Lyman (3 
decrement extending down to Lyman 7 at 972 A. The 7 and 5 decrements are also visible, 
although we will need more definitive statistics (below) to quantify their significance. The 
decrement below Lyman e at 937 A while present in the data, will turn out to be of 
questionable statistical significance in this sample. 

We can now quantify the observed relative absorptions of Lyman /3, 7, 6, and perhaps 
e, relative to Lyman a; in other words the ratios of the mean equivalent widths (the precise 
definition of which is equation [16] below) of different Lyman absorption lines in the clouds: 
Wp/W ai W 1 /W a , Ws/W a , etc. In the limited redshift range that we are studying, we will 
assume that the equivalent width ratios do not depend on redshift. These ratios are, as we 
will see in Paper II, direct indicators of the physical state of the clouds. 

An important point is that we are not simply quantifying the decrements seen in Figure 

5. That figure does not display the fact that each datum has an individual emission redshift 
z em associated with it. Use of these individual redshifts (through equation ||) allows, at least 
in principal, the removal of that part of the scatter in Figure 5 that is due to the range of 
redshifts in the quasar sample. For example, Wp/W a is estimated using all points with X em 
between 972 A and 1020 A. Each observed point generates a one-point estimate of the form 

Wp/Xp _ HQ- 1 ) -r a (\ obs /l2lQA-l) 
W a /X a ~ r a (A ofes /1025A - 1) 

where Q is the ratio of observed flux to the extrapolated continuum of the point (i.e., the 
observed transmission), and r = T(z ab& ) is given by equation (pj). Here, X a = 1216 A and 
A/3 = 1025 A are the laboratory wavelengths. The one-point estimates are then averaged. 

To estimate W 7 /W a , we use points in the range A em between 949 A and 972 A. Now 
we must subtract from the measured total optical depth ln(<5 _1 ) both Lyman a and Lyman 
(3 corrections at their respective absorption redshifts and relative strengths, 

W 7 /A 7 _ HQ- 1 ) - T a (X obs /121QA - 1) - r/3(A ofe /1025A - 1) 
W a /X a ~ r a {X obs /972A- 1) 
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where 

f a (z abs ) (14) 
is the derived mean optical depth in Lyman /3 as a function of redshift. 

One continues in this manner to obtain the ratio Ws/W a using points in the range 937 
A to 949 A, and W e /W a using points in the range 930 A to 937A. Obviously, looking at 
Figure 5, one must at some point begin to wonder whether the values obtained have any 
meaning. Also, it is clear that, because of the repeated subtraction of earlier ratios, errors 
in the determination of later ratios will be highly correlated (actually, anticorrelated) with 
the errors of earlier ratios. 

Once again, the bootstrap technique of resampling with replacement provides 
quantitative answers to these concerns. As for the previously determined parameters A and 
7, we have made 100 resampled determinations of the 4 ratios Wp/W ai W^/W a , Ws/W a , 
and W e /W a , choosing a different set of 29 quasars in each resampling. 

In the notation of the discussion leading to equation (|), let R i: i = 1, ... ,4, now 
denote the above four ratios, in the obvious order. Then the resulting measured quantities 
R®, a(Ri), Cij, Cy 1 , and are given in Table 3. Although the last ratio, W e /W a , is not 
well determined (detected at less than two standard deviations), the other ratios seem quite 
well established by the data. The first two ratios, Wp/W a and W^/W a , are determined by 
this data set with about 20% accuracy. We will see in Paper II that these values are able to 
constrain the nature of Lyman a clouds at high redshifts that are not yet accessible to high 
resolution studies. 



5. Fluctuations in Optical Depth around Mean Values 

In principle there is information not only in the mean values of Lyman a absorption 
(and the higher lines (3, 7, 5), but also in the statistics of the fluctuations of individual points 
around the mean. In Figure 5, for example, one wants to identify the mechanisms that 
contribute to the spread of the individual points around the mean, defined in that Figure 
to be unity between 1050 A and 1170 A. Measurement noise, which tells us nothing about 
the Lyman a clouds, is one contributing factor, as are uncertainties in our extrapolation of 
the quasar continua, which are also not of intrinsic interest. 

The more interesting contribution, because it potentially does give information about 
the clouds, is that of the Poisson statistics of how many clouds are present in each spectral 
resolution element. 
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Unfortunately, as we will now see, the exploitation of this information is not 
straightforward. 

Let us first revisit the question of the mean transmission. The absorption equivalent 
width W of a cloud (in its rest frame) depends on its physical parameters, for example 
column density N, temperature parameter b, and so on. Let p denote the vector of such 
parameters, and let jV(p)c?p be the number of clouds along the line of sight per unit rest 
wavelength in a volume dp of parameter space. Then the total number of clouds per unit 
rest wavelength is given by 

n = [ Af(p)dp (15) 



(if the integral converges), while the mean equivalent width of a cloud is 

W = - ! W(p)Af(p)dp (16) 
n J 

If the lines are randomly placed, e.g., without clustering, then the mean transmission Q is 
given by 

Q = exp(-nW) (17) 

if the integral for n exists, and 



Q = exp(- J W(pW(p)dp) (U 



otherwise (for example, if there are an infinite number of lines with negligibly small W). 
Equation (|i~7"f ) is so natural as to seem intuitively obvious, but it is in fact quite a nontrivial 
result, relating the average of a (nonlinear) exponential to the average of the exponent. A 
proof is given by Goody (1964, §4.5); see also the Appendices, equations (|A15| ), (|B24|) , and 

(HD- 

The sad fact is that there is no such universal result for the next moment, 

Var(Q)/g 2 = ( ( £ - 1 ) ) (19) 

2 

Rather, the variance Vax(Q)/Q depends in a nontrivial way on the joint distribution of 
equivalent widths (or, more fundamentally, column density N) and doppler widths b, as 
well as on the instrumental spectral resolution. 

Some simple model cases are derived in Appendix B and are illuminating. Assume a 
square line profile of full width Wq, and suppose that all lines have the same column density, 
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and thus the same equivalent width W. Let A be the instrumental spectral resolution. 
Then an exact result and its power series expansion are 

-2 



—2 2W n 
Var(Q)/Q 2 = 



A* 



where 



.4 



1 



nW 



nW 



A 

1+ 3 + 



(20) 



(21) 



See Appendix B, equations ( |B41| ) and ( |B44| ) . If the instrumental response is not square but 
is rather described by a response function w(X — Aq), then the definition of A is 



wdX 




w dX 



(22) 



(cf. equation [B30| in Appendix B). 



The limiting form of equation (|20|) for small A can be reproduced, up to a constant, by 
a simple physical argument: Because the lines have a width Wq, the instrumental width 
A contains A/Wo independent elements. In each of these, the mean optical depth is nW, 
while the number of lines contributing to this mean is tiWq. Thus, in each element, the 
r.m.s. fluctuation in optical depth, which is also the fractional fluctuation in transmission Q 
is about (nW)(nWo)~ l l 2 . Averaging over the independent elements reduces this fractional 
fluctuation by an additional factor (Wq/ A) 1//2 , yielding finally 



2 nW 

Var(Q)/Q 2 ~ — 



(23) 



The reason that the exact result ( P0|) has an additional exponentially increasing factor in A 
is that, for large optical depths and saturated lines, the fluctuations become dominated by 
the cases where increasingly rare "windows" happen to occur in the random placement of 
the lines. 



A second analytic result derived in Appendix B (cf. equations B41 and |B48| ) is for the 
case where there is no characteristic value for W, but rather a power-law distribution in 
N, the column density. (This is a more realistic idealization of the actual situation for the 
Lyman a clouds.) If the number of clouds along the line of sight with N between iV and 
iV + dN is proportional to N~@dN, then equation (p0|) is replaced by 



Var(Q)/g 2 2W ° 



A 



V 

K 



-1) 



where now 



k=(2- 2 p - 



ln= = (2 

Q 



2 /3- 



(24) 



(25) 
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(cf. equations |T7| and 18). When the mean optical depth r is not too large, we have 
approximately 

Var(g)/g 2 ^(2-2^ 1 )^r (26) 

A typical value for (3 is 1.5 (Sargent 1988a; we will have more to say about this value in 
Paper II). 

Let us now compare these results to the SSG data. We estimate Va.r(Q)/Q in the data 
by the following procedure, which is designed to greatly reduce the variance due to errors 
in the continuum fits: For each measurement Qi, we estimate Qi — Q not by subtracting 
the overall mean Q, but rather by subtracting a local weighted average of the Q^s for that 
particular quasar and neighboring wavelength bins. The weighted average is taken to have 
a triangular profile, with unit amplitude at bin i falling to zero at bins i ± 30 (that is, 
±300A in observed wavelength). With a mean redshift of about 1 + z = 4.5, this means 
that we are sensitive to variations in the extrapolated continuum only over scales of 65 A, 
which are likely negligible. Our measured result for the SSG sample is 

Var(Q)/Q 2 = 0.05 ± 0.01 (27) 

where the 1-cr error bar is again obtained by resampling. 

At a redshift 1 + z = 4.5 we have r = 0.67 (equation [| cf. Figure 4). Using SSG's 
quoted spectral resolution A = 25A, equation ([23]) implies Wo ~ 0.7A, or b ~ 80km s _1 . 
This value is implausibly large by about a factor of about 2. At this stage of analysis 
there are three possible resolutions: (1) Since we do not have an absolute calibration of 
the measurement noise for each data point, we cannot say with confidence that this excess 
variance is not an instrumental effect. (2) The excess variance might be due to a two-point 
correlation function £ in the clouds on velocity scales comparable to the instrumental 
resolution A, that is Ac/1216A ~ 6000km s^ 1 . Clustering of Lyman a clouds has previously 
been detected only at much smaller scales (e.g., Webb 1987, Crotts 1989). (3) The factor 
of 2 discrepancy might be an artifact of our simplistic analytic assumption of a square 
line profile; the observed fluctuations might in fact be consistent with a Poisson random 
distribution of clouds. 

In Paper II, further analysis will show that the third resolution is the most likely one. 
We note here that our finding fluctuations comparable to the predicted Poisson value helps 
confirm an implicit assumption that we have made throughout this paper: the overall 
absorption at high redshifts is not significantly due to a continuous Gunn-Peterson trough 
rather than a superposition of individual clouds. 



-16- 



6. Conclusions 

The principal results of this paper are a set of techniques for the statistical analysis 
of the Lyman a forest that can be used without the necessity of resolving and identifying 
individual lines. We have seen that there exist a set of well-defined parameters that can thus 
be measured with useful accuracy. Where there is overlap with previous results (including 
those obtained by high resolution studies) our methods give a reassuring agreement. 
Additionally, our methods make accessible some new parameters, such as the equivalent 
width ratios of the Lyman series. 

For the redshift range of about 2.5 < z a bs < 4.3 that is sampled by the lines of sight to 
29 SSG quasars, we have the following specific results: 

1. The mean Lyman a absorption at a redshift z is well approximated by 
T(z) = A(l + zY' +1 with 7 = 2.46 ± 0.37 and A = 0.0175 - 0.00567 ± 0.0002, in 
agreement with previous determinations at moderate redshifts from high spectral resolution, 
ground-based observations (though not with the low redshift measurements of HST). 

2. The ratio of mean Lyman (3 absorption to mean Lyman a absorption at a fixed 
absorption redshift (which is diagnostic of the physical state of the clouds at that redshift) 
is 0.476 ± 0.054. The ratios for Lyman 7, 5, and e are also measurable, and given, along 
with the error covariance matrix for all the ratios, in Table 3. 

3. Fluctuations around the mean absorption are comparable to, possibly a factor of 2 
larger than can be explained as, simple Poisson fluctuations in the number of clouds along 
the line of sight. The factor of 2 discrepancy could be due to measurement noise, model 
imprecision (see Paper II), or a small, positive two-point correlation function £ on a scale 
of 6000 km s _1 . The measurement of fluctuations comparable to the Poisson value argues 
against a significant Gunn-Peterson trough at high redshifts. 

We have benefited from discussions with John Bahcall and John Huchra. The referee 
made a number of helpful suggestions. This work was supported in part by the National 
Science Foundation (PHY-91-06678) and by NASA (NAG5-1618). 



A. Appendix A 

We here derive some statistical properties of the extinction optical depth function r(A) 
as a function of wavelength. We base our results on a simple statistical model, closely 
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related to random line models. Anticipating needs of future papers, our development is 
somewhat more general than strictly required for the present applications. 

We first express r(A) as a sum of contributions from a number of clouds along the line 
of sight, 



where T(A|Aj,pj) is the contribution from the ith cloud. This depends on the central 
wavelength of the cloud, Aj, as well as on a vector of parameters for the cloud, which 
may include its column density iVj, its velocity parameter fej, and possibly other parameters. 

It may be helpful at this point to remark that the particular form for T used in the 
body of the paper is T(A|Aj,pj) = iVja:(A — Aj,fej), where iVj is the column density, 6j the 
velocity parameter, and a(X — Aj, bi) is the atomic absorption coefficient, typically a Voigt 
function. In this special case A, simply defines the central wavelength of the profile, and the 
parameters pj consist only of JVj and 6j. We shall generally assume that the A dependence 
of T falls off sufficiently rapidly away from Aj to assure the convergence of certain integrals. 

It is easiest to express our statistical model as a two-step process: First of all, we 
assume that the positions of the central wavelengths Aj of the clouds are Poisson random 
distributed with a given mean density of lines per unit wavelength interval. Second, given 
the positions Ai, A 2 , . . ., of the clouds, the parameters pi, p 2 , . . . are assumed to be 
independently random, that is, their joint conditional distribution function is a product of 
independent factors, one for each cloud, 



Thus the distribution of the parameters pj of the ith cloud can depend on the associated 
Aj, but not on the p's or A's of any other clouds. 

Let the mean density of line centers per unit wavelength interval and in a volume 
dp of parameter space be given by A/"(A, p)dp. Generalizing the discussion in §5, we have 
here allowed the line density function M to depend on wavelength. For convenience in the 
following derivations we shall assume that the mean total density of lines at A, 



is finite, temporarily ignoring the possibility (alluded to in the text) of a power law 
divergence at small column densities. We may then take 



r(A)=^T(A|A,, Pi ), 



(Al) 




(A2) 




(A3) 




(A4) 
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as the conditional probability density for the parameters p, for the zth cloud. 

The statistical model defined above is similar to the well-known Poisson or shot process 
(see, e.g., Parzen 1962). However, this process has been generalized here to allow the 
mean density of centers to vary with A and the function T to depend on A and \ in more 
complicated ways than simply through their difference. This generalization is not strongly 
needed for the applications of this paper, but will play a more important role for later 
papers in this series. 

A very powerful and general way to express the statistical properties of the function 
r(A) is through its characteristic functional (see, e.g., Bartlett 1955). This is defined for an 
arbitrary function //(A) by 



= (exp 



i / / u(A)r(A) dX 



(A5) 



where the angle brackets denote averaging over the distribution of line centers and over the 
parameter space p for each line. 



The range of integration in equation ( pAq ) has been left unspecified. In principle we 
would like it to be the infinite range of A, but for purposes of the following derivation it is 
convenient to choose it temporarily to be some definite, large (but finite) interval, which 
adequately spans the essential nonzero region of the function //(A). Such a finite region of 
integration implies that that the mean number of clouds in the interval 

M = J n{X)d\ = J dxj dpA/"(A,p). (A6) 

is finite. The actual number of lines M in the interval is a random variable with Poisson 
distribution function 

Pu = ^j-e-, (A7) 



which can range from to oo. For each realization of M, the summation in equation flAy ) 
is from i = 1 to M. 

We may now evaluate the part of the averaging in equation ( |A5| ) due to the Poisson 
distribution of line centers. We note that the probability that a line occurs in a differential 
interval dXi about Aj is simply dXin(Xi)/M, and that the probability distribution for each 
line is independent of all other lines (this is the Poisson assumption). Then 




M=0 



M 

I 

8=1 



M r 

H J MA)T(A|Ai,pi)dA 



(Ah 
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where now the angular brackets refer only to the average over all the parameters pi, P2, 
. . .. With equations ( |A2| ) and (^), we have, 



— 00 M 



£ -IK III [ dXi ~W~ exp \ l / MA)T(A|A,, Pi) d\ 

W=0 ' i=l \ m L 

z / /i(A)T(A|Aj, p^ dA 



■ exp 



(A9) 



We note that all factors in the product are identical, since the dummy variable of integration 
is irrelevant, so the product reduces to a simple power of one factor, which we may choose 
to be the one for % = 1, 



1 



M=0 



$[//] = e- M £ — { I d\ 1 I dp^A^p^exp 



i / /i(A)T(A|Ai, pi) (iA 



(A10) 



Then using equation ( |A^ ) and the power series for the exponential, along with some changes 
in dummy variables of integration, we find 



exp 



i / MAX A') dX 1 



expj— J dX J dpAf(X,p) ^1 — exp % J fi(X')T(X'\X, p) dX' 



(All) 



Recall that finite limits of integration were introduced to make the quantity M finite. 
However, assuming a sufficiently localized /i(A), such divergent quantities no longer appear 
in equation ( All ), so the integrals in this formula can be considered to have arbitrary 
limits. Similarly, we note that the quantity ^(A) no longer appears, so in many cases this 
equation can have meaning even when the integral in equation (|A3|) diverges. 



Equation ( |A11| ) is the major result describing the statistical properties of the underlying 
Poisson process. We can use it to derive a variety of important statistical results. For 
example, if the two sides of the equation are expanded to first order in fi(X') and the 
corresponding terms equated, one obtains the result for the mean of t(A'), 



(t(A')> = Jdxj dpAf(X,p)T(X'\,X,p). 



(A12) 



Equation ( |A11|) can also be expanded to second order. In doing so, one must be careful 
to introduce a new dummy variable of integration A" on one factor. Omitting details, we 
obtain the correlation function, 

([r(A') - (r(A'))] [r(A") - (r(A")>]> = / dX J dp M(X, p)T(X'\, X, p)T(X"\, X, p). (A13) 
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This result generalizes the well-known Campbell's theorem (see, e.g., Parzen 1962) by 
allowing A dependence in the statistical quantities. 

For the study of Lyman a clouds, more important than averages involving r itself are 
averages involving the exponential extinction law, 



q(\) = e- T(A) . 

Setting /x(A') = i8(X' — A) in equation ( |A11|) gives immediately, 



q(X) = (q(X)) = exp - / d\' / dpj\f(\, p 



-T(A'|A,p) 



(A14) 



(A15) 



This formula represents a generalization of equation ( |IBD of the text. 

Correlation properties of q(X) can be found similarly. Setting //(A A 
iS(X"' - A') + i5{X"' - X") in equation (|ATTD gives 



(q(X')q(X")} = exp [- J dX" 1 J dpJ\f(X M , p) (l - fWWM) 



(A16) 



[We note that a special case of equation ( |A16|) , for A' = A", and for a homogeneous model, 
was given by M0ller, P., and Jakobsen (1990)]. Using equation ( |A15| ) twice, this can be 
written 

(A17) 



(q(X')q(X")) = q(X')q(X")e 



where 



H(X', A") = J dX'" J dpM{X"', p) (l - e - T ( A 'l A "' p )) (l - e - T ( A "l A "' p )) . (A18) 



We may also write 

([q(X') - q(X')\ [q(X") - q(X")}} = q(X')q(X") (e^( A '' A ") - l) . 



(A19) 



Note that H(X', A") is small when the wavelength difference |A' — A"| is much larger than 
the atomic line widths, since then the T functions in equation ( |A18| ) do not significantly 
overlap. Thus equation ( |A19| ) shows that the values of q(X) are essentially uncorrelated for 
such wavelength differences. 



B. Appendix B 
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The purpose of this Appendix is to derive equations ( p0|) through (p5|) of §5, applying 
the general results of Appendix A. 

Before proceeding, we want to take into account that the measured extinction is not 
q(X), but rather an average over the instrumental profile, 

Q(A) = fw{\- A')g(A') dX'. (B20) 

It is usually sufficiently accurate to take w to be a function of wavelength differences A — A'. 
In this Appendix, though not in the main text, it is defined with normalization 

1 = f w (X-X')dX'. (B21) 



Besides Q(X) itself, we are also interested in its square 

Q 2 (A) = / dX'w(X - A') / dX"w(X - X")q(X')q(X"), (B22) 



which is needed in defining the variance of Q. A useful generalization of this is the product 
of Q's at two different points, 

Q(Xi)Q(X 2 ) = J dX'w(X l - A') J dX"w(X 2 - X")q(X')q(X"), (B23) 

which is needed in defining the correlation properties of Q. (The variables Ai and A2 here 
are general variables, and have nothing to do with the i = 1 and % = 2 clouds.) 

We shall now assume that the distribution Af is independent of A (or at least that its 
scale of variation in A is much larger than either the atomic line profile or the instrumental 
resolution). We also assume that T(A|Aj,pj) = A^a(A — Aj,6j), which depends on A and Aj 
only through their difference A — Aj, at least locally. 



Substituting this expression for T into equation (|A15|) , we immediately obtain, 



where 



(g(A )> = (e^) = exp (- J W(p)M(p) dp) , (B24) 
W(p) = J dX {1 - exp [-Na(\, b)}} . (B25) 



is the equivalent width of a line with column density iV and velocity parameter b. In 
deriving equation ( |A13|) , we have shifted the variable of integration, showing that VK(p) 
is independent of Ao (at least over wavelength intervals sufficiently small that a can be 
considered to be only a function of wavelength differences). Consequently (q) is also 
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independent of A . Then from equations (|B20|) and (|B21| ), it follows that Q = (Q) = (q). 
Thus equations ( |B24| ) and ( |B25| ) constitute an independent proof of equation ( |T8| ) of the 
text. 

We next consider the problem of determining the statistical average of Q 2 , given in 
equation (|B22|) . Actually, it is almost as easy to determine the average of the product 
of <5's in equation ( |B23| ), so we shall do this. First of all, we need the average of 
exp[— t(A') — t(A")]. This is given by equation (|A19| ), which can now be written 



([g(A')-g][g(A")-g]) = g 2 ( e - H(A "- A,) -i). 



(B26) 



where 



H(X) = J dpM(p) J d\'" (l - e - N ^ x '"^) (l - e - Na ^'" +x ^) . (B27) 

Here we have made appropriate changes of variables to take advantage of the difference 
dependence of the atomic absoption coefficients. Using equation ( |B23| ), we find that 

([Q(Xi)-Q] [Q(X 2 )-Q}) = Q 2 J dX'w(X 1 — A') J dX" w(X 2 - A' + A") {e~ H ^ - l) . 

(B28) 

Equation ( |B28| ) can be substantially simplified under conditions where the individual 
atomic lines are very much under-resolved, that is, when the instrumental function W is 
much broader than atomic widths. (This is true for the data of the present work.) In this 
case, the effective range of the A" integration is very much smaller than the width of the 
instrumental profile, and we can ignore the A" dependence in the function w(X\ — X' + A"), 
setting it equal to u>(Ai — A') and taking it outside the A" integration. This gives 



Q(Ai)-Q Q(A 2 )-Q 



Q 2 



dX" w(X")w(X 2 - X± + X") 



dX" 



-H{\") 



- 1 



(B29) 



The quantity in brackets determines the dependence of the correlation function on the 
wavelength difference Ai — A. When Ai = A, the LHS of this equation is equal to the 
variance of Q, denoted Var(Q). Thus 



Var(Q) = Q 2 



f(A") 



dX"w\X") JdX" (e- H 
Thus equation ( |B29|) can be written in the simple form 

<5(Ai) - Q] \Q(X 2 ) -Q})= Var(g)F (A 2 - Ai) 



where 



Y(X) 



JdX"w(X")w(X + X") 
JdX"w 2 (X") 



(B30) 

(B31) 
(B32) 
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For the data of this paper, the instrumental function is well represented as the 
rectangular function, 

u ,(A) = (f'' (B33) 
v ; \ 0, A > A/2. v ; 

In this case, 

Y(X) = (1 - |A|/A) + (B34) 

where (•••)+ implies zero for negative arguments. This shows that the correlation function 
has a central value equal to the variance, and has a triangular shape, going to zero for lags 
greater than the instrumental width. 

It remains only to determine the variance Var(Q) itself. For the rectangular 
instrumental profile this is given by 

Var (Q)/Q 2 = ^J dX" (e~ H ^ - l) . (B35) 

In general the integral in equation ( |B35|) must be evaluated numerically, which also 
involves numerical evaluation of the function H(\) defined in equation (|B27|). However, 
considerable insight can be obtained by invesigating the simple model using the single 
rectangular line profile, 



a(A-A',6)= « " ? ' (B36) 



f OtQ, 


|A 


-A'| 


1°. 




-A'| 



where «o is the height of the profile and Wq is its full-width. The parameter b does not 
appear here, since it implicitly takes a single value, which sets the values of the constants 
a>o and Wq. The parameter vector p then consists only of the column density N. 

With this choice of profile equation ( |B25| ) gives 

W(N) = Wo (l - e"^ ) , (B37) 

and equation (|B24j) then gives 



Q = exp 



-W J dNAf(N) (l-e 



-Nat 



(B38) 



Having found the mean transmission, we next want to find the variance (and correlation 
function). From equation ( |B27| ) it follows that, 

F(A) = k(1-|A|/W ) + , (B39) 
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where 



-Nag 



k = W J dNAf(N) (l-e 
The variance is then found from equation (P35| ), 

var(g)/g 2 = 1 r dx ( e -«a-Wo) _ i) 2W ° 

A Jo v ' 



A 



- (e K - 1) - 1 

K 



(B40) 



(B41) 



The simplest special case of these results is when there is also a single column density 
N for all clouds. Is this case J\f(N) = nS(N — No), so that equations ( |B37| ) and ( |B38| ) 
imply that the mean transmission is 



Q = exp —nW 
where the equivalent width of each line is, 

W = W (l - e 

The variance is given by (|B41 ), where now, 



-N a 



(B42) 



(B43) 



k = nW (l-e 



-Na 



nW 
~Wo~ 



(B44) 



which proves equation (21) of the text. 



We can also find results for the case where Af is the power law in N, J\f(N) = KN 13 . 



Then 



where 



Q = exp 



-nW 



(B45) 



W 



KWqJ dNN- p (l-e 



-Na 



KW a 



dN N 1 ~ l3 e~ Na ° 



(3-1 Jo 
= KW at 1 T(2-(3)/((3-l) 

An integration by parts has been used, along with the definition of the T function. 

Similarly, we now find from equation ( |B40| ). 

-Na 



(B46) 



POO . 

k = nKWo / dN N~P (l-e 
2nKW °ao f" ^ ( 



(l Na _ e -2Afao 



(3-1 Jo 
(2 - 2l i - l )nKWoat l T{2 - (3)/((3 - 1) 



(B47) 
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Using equation ( |B46|) we now find 



k= (2-2 1 ~ l3 )nW = (2-2 /3 - 1 )ln^, (B48) 



which proves equation (|25j ) of the text 
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Line 


A (A) 


w (A) 


Lja 


1216. 


20. 


NV 


1240. 


18. 


01 


1302. 


7.2 


CII 


1335. 


6.0 


SilV/OIV 


1400. 


11. 


CIV 


1549. 


15. 


Hell/OIII 


1651. 


77. 


Fell feature 


1770. 


26. 


A1III 


1858. 


30. 


cm 


1909. 


22. 



Table 2: Lines, emission wavelengths, and nominal widths used for fitting the underlying 
continua of individual quasars in the range 1250 A to 2200 A. Nominal widths are used only 
as the starting point for a linearized width correction. 
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Ratio 


w p /w a 


WJW a 


W s /W a 


w e /w a 


i = 


1 


2 


3 


4 


Measured Value 


0.476 


0.351 


0.157 


0.118 


Standard Error ( 1-a ) 


0.054 


0.072 


0.082 


0.072 


Correlation Matrix Cij 










3 = 1 


0.00296 


-0.00156 


-0.00005 


0.00085 


J = 2 


-0.00156 


0.00532 


-0.00266 


-0.00041 


J = 3 


-0.00005 


-0.00266 


0.00668 


-0.00351 


J =4 


0.00085 


-0.00041 


-0.00351 


0.00522 


Inverse Correlation Matrix C,-,- 1 








3 = 1 


422. 


147. 


49.5 


-23.6 


J = 2 


147. 


355. 


224. 


155. 


J = 3 


49.5 


224. 


378. 


264. 


J =4 


-23.6 


155. 


264. 


385. 


Coefficient of Correlation 








3 = 1 


1.00 


-0.39 


-0.01 


0.22 


3 = 2 


-0.39 


1.00 


-0.45 


-0.08 


J = 3 


-0.01 


-0.45 


1.00 


-0.59 


J =4 


0.22 


-0.08 


-0.59 


1.00 



Table 3: Measured quantities associated with the ratios of mean equivalent widths for Lyman 
a, (3, 7, 5, and e in the SSG sample. In the column headings W a = W a /\ a , Wp = Wp/\p, 
etc. 



-30- 



Fig. 1. — The observed spectra of 29 SSG quasars are here superposed after shifting each 
to its emission rest frame and scaling each to a common magnitude at 1450 A. Despite 
indisputable differences in the individual quasars' continuum slopes and emission features, 
there is considerable similarity in the spectra. The principal interest of this paper is in the 
statistical analysis of the Lyman a forest shortward of 1200 A. 

Fig. 2. — Fitted line profiles and derived continuum models for each of the 29 SSG quasars 
analyzed. Fitting is done by a linear model, which can give artifacts in line wings, but is 
otherwise more robust than a general nonlinear fit (see text). The purpose of these fits is to 
obtain continuum models that can be extrapolated shortward of 1200 A. 

Fig. 3. — Extrapolations of the continuum models shown in Figure 2 to emitted wavelengths 
between 930 A and 1200 A. Each point shown, taken as a fraction of the extrapolated 
continuum above it, is an (approximately independent) measurement of Lyman a absorption 
at a calculable redshift. The ensemble of all the points in this Figure (excluding a small range 
of emitted wavelengths around Lyman (3) is the data set that is analysed in the rest of this 
paper. 

Fig. 4. — Left, the subset of points from Figure 3 with emitted wavelength between 1050 A 
and 1170 A (each normalized to its extrapolated continuum level) is plotted as a function of 
emitted wavelength. Since the 29 SSG quasars vary substantially in redshift, the observed 
transmission varies widely, with no significant trend. Right, the same data is plotted as 
a function of observed wavelength or, equivalently, absorption redshift for Lyman a (top 
scale) . Here the trend with redshift is clear. The solid line fits a power law model with mean 
optical depth varying as (1 + z) 1+1 ) with 7 = 1.46. The shaded band approximates the range 
of statistical uncertainty in the fit, as determined by the bootstrap method of resampling 
the 29 quasars (see text). 

Fig. 5. — The individual transmission measurements are shown as a function of emitted 
wavelength, after correcting each point for Lyman a absorption at its own absorption 
redshift. The shaded band is a moving average of the points. One sees the Ly-/3, Ly-7, 
Ly-<5, and possibly Ly-e decrements. These are jointly fitted, after shifting each point in the 
proper emission region to its proper absorption redshift, and equivalent width ratios for the 
corresponding Lyman lines, along with an error covariance matrix (again obtained by the 
bootstrap method), are obtained. See text for details. These mean equivalent width ratios 
place significant constraints on physical conditions in the clouds. 



